install.packages("PanelMatch")
library(PanelMatch)
rm(list=ls())

## Set the working directory to the folder where the dataset file "visit_defspending-simple.csv" is located
#setwd("/Users/jamesdkim/Library/CloudStorage/Dropbox/Working Papers/Leadership Visits and Military Expenditures/FPA Second R&R/Replication files submission")

## Import the dataset file
#visit_defspending.simple <- read.csv("~/Library/CloudStorage/Dropbox/Working Papers/Leadership Visits and Military Expenditures/FPA Second R&R/Replication files submission/visit_defspending-simple.csv")

# get a matched set without refinement
sets0 <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                    treatment = "visit_nolastyear_official", refinement.method = "none",
                    data = visit_defspending.simple, match.missing = FALSE,
                    covs.formula = ~ visit_nolastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                    size.match = 3, qoi = "att",
                    outcome.var = "defspend_ch",
                    lead = 0:4, forbid.treatment.reversal = FALSE)

# get another matched set with refinement using ps.weight
sets2 <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                    treatment = "visit_nolastyear_official", refinement.method = "ps.weight",
                    data = visit_defspending.simple, match.missing = FALSE,
                    covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                    size.match = 3, qoi = "att",
                    outcome.var = "defspend_ch",
                    lead = 0:4, forbid.treatment.reversal = FALSE)

# get another matched set with refinement using ps.match
sets3 <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                    treatment = "visit_nolastyear_official", refinement.method = "ps.match",
                    data = visit_defspending.simple, match.missing = FALSE,
                    covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                    size.match = 3, qoi = "att",
                    outcome.var = "defspend_ch",
                    lead = 0:4, forbid.treatment.reversal = FALSE)

# get another matched set with refinement using CBPS.weight
sets4 <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                    treatment = "visit_nolastyear_official", refinement.method = "CBPS.weight",
                    data = visit_defspending.simple, match.missing = FALSE,
                    covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                    size.match = 3, qoi = "att",
                    outcome.var = "defspend_ch",
                    lead = 0:4, forbid.treatment.reversal = FALSE)

# get another matched set with refinement using CBPS.match
sets5 <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                    treatment = "visit_nolastyear_official", refinement.method = "CBPS.match",
                    data = visit_defspending.simple, match.missing = FALSE,
                    covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                    size.match = 3, qoi = "att",
                    outcome.var = "defspend_ch",
                    lead = 0:4, forbid.treatment.reversal = FALSE)

# Figure 2 in the online appendix
balance_scatter(non_refined_set = sets0$att,
                matched_set_list = list(sets2$att, sets3$att),
                data = visit_defspending.simple,
                covariates = c("visit_lameduck_official", "realgdp", "ln_cinc", "demo", "midcount", "l_postcw", "landcontig", "viol", "rivalry", "growth"))

legend(x = 0, y = 0.8,
       legend = c("PS weighting",
                  "PS matching"),
       y.intersp = 0.65,
       x.intersp = 0.3,
       xjust = 0,
       pch = c(1, 3), pt.cex = 1,
       bty = "n", ncol = 1, cex = 1, bg = "white")


# Figure 3 in the online appendix
balance_scatter(non_refined_set = sets0$att,
                matched_set_list = list(sets4$att, sets5$att),
                data = visit_defspending.simple,
                covariates = c("visit_lameduck_official", "realgdp", "ln_cinc", "demo", "midcount", "l_postcw", "landcontig", "viol", "rivalry", "growth"))

legend(x = 0, y = 0.8,
       legend = c("CBPS weighting",
                  "CBPS matching"),
       y.intersp = 0.65,
       x.intersp = 0.3,
       xjust = 0,
       pch = c(1, 3), pt.cex = 1,
       bty = "n", ncol = 1, cex = 1, bg = "white")


# Figure 4 in the appendix: Estimated Effects of Treatment Over Time
PM.results <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                         treatment = "visit_nolastyear_official", refinement.method = "ps.match",
                         data = visit_defspending.simple, match.missing = TRUE,
                         covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                         size.match = 5, qoi = "att",
                         outcome.var = "defspend_ch", lead = 0:5, forbid.treatment.reversal = FALSE)

PE.results <- PanelEstimate(sets = PM.results, data = visit_defspending.simple)
plot(PE.results)


# Figure 5 in the appendix: Estimated Effects of Treatment Over Time
PM.results <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                         treatment = "visit_nolastyear_official", refinement.method = "ps.match",
                         data = visit_defspending.simple, match.missing = TRUE,
                         covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                         size.match = 5, qoi = "att",
                         outcome.var = "defspend_ch", lead = 0:5, forbid.treatment.reversal = FALSE)

PE.results <- PanelEstimate(sets = PM.results, data = visit_defspending.simple)
plot(PE.results)


# Figure 6 in the appendix: Estimated Effects of Treatment Over Time
PM.results <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                         treatment = "visit_nolastyear_official", refinement.method = "CBPS.match",
                         data = visit_defspending.simple, match.missing = TRUE,
                         covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                         size.match = 5, qoi = "att",
                         outcome.var = "defspend_ch", lead = 0:5, forbid.treatment.reversal = FALSE)

PE.results <- PanelEstimate(sets = PM.results, data = visit_defspending.simple)
plot(PE.results)

# Figure 7 in the appendix: Estimated Effects of Treatment Over Time
PM.results <- PanelMatch(lag = 2, time.id = "year", unit.id = "ccode",
                         treatment = "visit_nolastyear_official", refinement.method = "CBPS.weight",
                         data = visit_defspending.simple, match.missing = TRUE,
                         covs.formula = ~ visit_lastyear_official + realgdp + ln_cinc + demo + I(lag(midcount, 1)) + l_postcw + landcontig + viol + rivalry +growth,
                         size.match = 5, qoi = "att",
                         outcome.var = "defspend_ch", lead = 0:5, forbid.treatment.reversal = FALSE)

PE.results <- PanelEstimate(sets = PM.results, data = visit_defspending.simple)
plot(PE.results)
